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Quantum ground-state problems are computationally hard problems; for general many-body 
Hamiltonians, there is no classical or quantum algorithm known to be able to solve them efficiently. 
Nevertheless, if a trial wavefunction approximating the ground state is available, as often happens 
for many problems in physics and chemistry, a quantum computer could employ this trial wavefunc- 
tion to project the ground state by means of the phase estimation algorithm (PEA). We performed 
an experimental realization of this idea by implementing a variational-wavefunction approach to 
solve the ground-state problem of the Heisenberg spin model with an NMR quantum simulator. 
Our iterative phase estimation procedure yields a high accuracy for the eigenenergies (to the 10~^ 
decimal digit). The ground-state fidelity was distilled to be more than 80%, and the singlet-to-triplet 
switching near the critical field is reliably captured. This result shows that quantum simulators can 
better leverage classical trial wave functions than classical computers. 

PACS numbers: 03.67.Lx 07.57.Pt 42.50.Dv 76.60.-k 



I. INTRODUCTION 

Quantum computers can solve many problems much 
more efficiently than a classical computer [1 . One gen- 
eral class of such problems is known as quantum simu- 
lation [H [3]. In this class of algorithms, the quantum 
states of physical interest are represented by the quan- 
tum state of a register of controllable qubits (or qudits), 
which contains the quantum information of the simulated 
system. In particular, one of the most challenging prob- 
lems in quantum simulation is the ground-state prepara- 
tion problem [4 of certain Hamiltonians, which can 
be either classical or quantum mechanical. Remarkably, 
every quantum circuit [5j , and even thermal states [H |7] , 
can be encoded into the ground state of certain Hamilto- 
nians, and purely mathematical problems, such as factor- 
ing [8] , can also be solved by a mapping to a ground-state 
problem. 

On the other hand, the ground-state problem has pro- 
found implications in the theory of computational com- 
plexity [9]. For example, finding the ground-state of a 
general classical Hamiltonian (e.g. the Ising model) is in 
the class of N P (nondeterministic polynomial time) com- 
putational problems, meaning that while finding the so- 
lution may be difficult, but verifying it is efficient when 
employing a classical computer. The Ising model with 
non-uniform couplings is an example of an NP-problem 
(more precisely, NP-complete) [10 . The quantum gener- 
alization of NP is called QMA (Quantum Merlin Arthur) 
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[5l. In this class, the verification process requires a quan- 
tum computer, instead of a classical computer. An ex- 
ample of a problem in QMA is the determination of the 
ground-state energy of quantum Hamiltonians with two- 
body (or more) interaction terms [11]. So far, there is no 
known algorithm, classical or quantum, that can solve all 
problems efficiently in NP and QMA. 

Most of the problems in physics and chemistry, how- 
ever, exhibit special structures and symmetries, that 
leads to methods for approximating the ground state with 
trial states |?/^t) (or trial wave-function) possible. For 
example, in quantum chemistry the Hartree-Fock 

mean field solution often captures the essential informa- 
tion of the ground state |eo) for a wide range of molecular 
structures. However, the applicability of these trial states 
will break down whenever the fidelity, 

l(eo|^r)|' , (1) 

quantified by the square of the overlap between the trial 
state \iI^t) and the exact state |eo), is vanishingly small. 
Specifically, if the fidelity of a certain trial state for a par- 
ticular many-body problem is small, for example, about 
F = 0.01, it might be considered as a "poor" approxi- 
mation to the exact ground state [13], when used as an 
input state in classical computation. For quantum com- 
puting, however, the same trial state can be a "good" 
input, as one only needs to repeat the ground-state pro- 
jection algorithm, e.g., by Abrams and Lloyd [M] (see 
below), for about O(IOO) times, which is computationally 
efficient especially when the Hilbert space of the many- 
body Hamiltonian is usually exponentially large. This is 
the motivation behind our experimental work. 

Several theoretical studies [15 -18 along this line of rea- 
soning have been carried out for various molecular struc- 
tures. Here we performed an experimental realization of 
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this idea with one of the simplest, yet non-trivial, phys- 
ical systems, namely the Heisenberg spin model in an 
external field. Our goals for this study are: (i) to de- 
termine the eigenvalues of the ground state, and (ii) to 
maximize, or to distill, as much as possible the ground- 
state from a trial state, which contains a finite {F = 0.5) 
ground-state fidelity. For (i), we employed a revised ver- 
sion of the iterative phase estimation procedure to de- 
termine the eigenvalues of the Hamiltonian (to the 10~^ 
decimal digit). Subsequently, we apply a state-filtering 
method to extract the ground-state fidelity from the final 
state to achieved (ii). For this study, we specifically chose 
three cases corresponding to three different values of ex- 
ternal field in the simulation, namely h = 0^ h = O.TS/ic, 
and h = 1.25hc^ where he is the critical value of the exter- 
nal field at which the ground-state and the first excited 
state cross each other (see Fig.[T]). This is a singlet-triplet 
switching, and our experimental simulation captures the 
change of the ground state around this critical point re- 
liably. 

Finally, we note that the approach employed here is dif- 
ferent from the method for preparing many-body ground 
states based on the adiabatic evolution p^H26] , where the 
initial state is usually chosen as the ground state of some 
simple Hamiltonian, which can be prepared efficiently, 
instead of the trial states, which aim to capture the es- 
sential physics of the exact ground state. The perfor- 
mance (complexity) of the adiabatic approach depends 
on the energy gap along the entire evolution path. In 
our approach, the performance depends on the fidelity of 
the initial state and the energy gap of the Hamiltonian. 
Furthermore, in these experiments (except Ref. [23 ), the 
eigen-energy and the ground state of the Hamiltonian 
are not usually determined simultaneously, and therefore, 
cannot be considered as completely solving the ground- 
state problem [4]. In spite of the differences between 
these two approaches, it is possible that the adiabatic 
method can be incorporated in our procedure to further 
enhance the ground-state fidelity of the final state. How- 
ever, this possible extension is not considered here. 

This paper is organized as follows: first we will provide 
the theoretical background for this experimental work. 
Then, we define the Hamiltonian to be simulated and the 
choice and the optimization of the initial state. Next, 
we outline the experimental procedures, and explain a 
revised iterative phase estimation algorithm. Finally, the 
experimental results will be presented and analyzed by 
a full quantum state tomography. We conclude with a 
discussion of the results and the sources of errors. 



II. THEORETICAL BACKGROUND 

The central idea behind this experimental work has 
a counterpart in the time-domain classical simulation 
methods [27 . In the context of quantum computing, 
the method was introduced by Abrams and Lloyd jlj. 
Specifically, it was shown that for any quantum state 
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FIG. 1: (Color online) (a) The energy eigenvalues versus ex- 
ternal magnetic field of the Heisenberg Hamiltonian (defined 
in Eq. The optimized state |V^*) = IV^ (-7r/4, 7r/2)) 

(see Eq. T4|) (black line) contains a linear combinations of 
the two eigenstates (red and blue lines) only, (b) The three- 
qubit NMR quantum simulator consists of a sample of ^Re- 
labeled Diethyl-fluoromalonate dissolved in ^H- labeled chlo- 
roform. The nuclear spins (circled) of ^^C and ^H are used 
as the system qubits and that of ^^F is the probe qubit. The 
parameters of the NMR couplings of this molecule are listed 
in the table. 



IV^) — ^k^kl^k) which has a finite overlaps \ak\ (or 
fidelity) with the eigenstates \ek) of a simulated Hamil- 
tonian, the phase estimation algorithm [2F will map, 
with high probability, the corresponding eigenvalues to 
the states of an ancilla quantum register. 



(2) 



Consequently, a projective measurement on the register 
qubits will, ideally, collapse the quantum state of the sys- 
tem qubits into one of the eigenstates. By analyzing the 
measurement outcome, one can determine the ground- 
state eigenvalue Eq, and even project the exact ground 
state |eo). 

Given any trial state IV^t), the performance of the al- 
gorithm depends on the overlap |aop, which can be max- 
imized using many classical methods, such as using ad- 
vanced basis sets matrix product states (MPS) rep- 
resentations [30], or any suitable variational method. 
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FIG. 2: (Color online) (a) The quantum circuit diagram for the experiment. The explicit construction of the unitary operators 
W and V{t) is detailed in the Appendix [31 . The boxed quantum gates in (a) (and pulse sequences in (b)) generate the input 
state which is an optimized variational state with respect to the Hamiltonian defined in Eq. (|3|. (b) The entire pulse 

sequence corresponds to the quantum circuit diagram for the case of zero external field, /i = 0. The complexity and the lengths 
of the pulse sequences for the other cases, namely h = 0.75/ic and h = 1.25/ic, are roughly the same as this pulse sequence. 



III. THE HAMILTONIAN AND THE 
OPTIMIZED INPUT STATE 

The method proposed here can be generalized to apply 
to more general Hamiltonians, but as a very good exam- 
ple, we will employ the Heisenberg Hamiltonian with an 
external magnetic field pointing along the 2;-direction: 



H = j {i^jl + ipl + i^il) + h {It + /: 



(3) 



where ^cr^, and cr^ is one of the Pauli matrices 

(a = X, z) acting on the k = a^b spin. On the other 
hand, in general, there is no restriction to the choice of a 
trial state, as long as it is not orthogonal to the ground 
state (in this case, the ground state algorithm necessarily 
fails) . To mimic the behavior of the commonly-employed 
trial states of more general systems, we require our trial 
state to satisfy the following conditions: (a) that it con- 
tains one or more parameters which can be adjusted to 
minimize the energy (H) , and that this procedure usually 
does not lead to the exact ground state, and (b) that it 
may capture only part of the vector space spanned by the 
eigenstates of the Hamiltonian H. One possible choice 
that fulfills the above criteria is the following variational 
state which contains two adjustable parameters, 6 and 0, 
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Here, \0) = cos (9 1 10) +sinl9|01) and \(p) = cos(p\00) + 
sin(^|ll). In general, the optimized states for each 
given pair of (J, h) are not necessarily the same. How- 
ever, in our case, we found that the optimized state 
Itp^) = 1?/^ (— 7r/4, 7r/2)) can minimize the energy for all 
values of h and J > 0. Moreover, it turns out that 
this optimized state captured two out of the four eigen- 
energies (see Fig.^) only; therefore, a single probe qubit 



is sufficient to resolve them (for more general cases, see 
the Appendix [31 ). We note that the fidelity F (cf. 
Eq. ([T])) of the state |?/^*) with the exact ground state 
|eo) is exactly 50%. 



IV. OUTLINE OF THE METHOD 

This algorithm starts with a set of system qubits ini- 
tialized in the state = and a single 
"probe" qubit in the (|0) + |1)) /V2 state. For differ- 
ent times t, a controlled U (t) gate, where U (t) = e~^^^ 
{h = 1), is then applied, resulting in the following 



state: (l/V^) (|0) + e''^^' |1)) |efc), where uj^ 

Ek. The reduced density matrix of the probe qubit, 

contains the information about the eigenvalues in its off- 
diagonal matrix elements, which can be measured effi- 
ciently in an NMR setup (see Appendix ^J). A classical 
Fourier analysis on the off-diagonal matrix elements at 
different times yields both the eigenvalues Uk and the 
overlaps \akf . To obtain the value of cj^ with high 
accuracy, a long time evolution of the simulated quan- 
tum state is usually needed. However, for Hamiltonians 
with certain symmetries, we can perform a simplified ver- 
sion of the iterative phase estimation algorithm (IPEA), 
which is similar but not identical to the ones performed 
previously in Ref. [23l |32] . We will explain the details of 
this IPEA in Section ED 

Once the ground state eigenvalue Eq of the Hamilto- 
nian H is determined, one can, for example, employ the 
state-filtering method [33| to isolate the corresponding 
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FIG. 3: (Color online) The absolute amplitude of the eigen- 
value spectra g(E) of the Heisenberg Hamiltonian as defined 
m Eq. These are obtained for three different values of 

the external magnetic field h, namely h — 0^ h — 0.75hc, and 
h = 1.25/ic, where he is the critical field at which the ground 
state becomes degenerate. The shaded region highlights the 
location the the singlet state depicted in Fig. [l] The arrows 
indicate the direction of the peak shift when the simulated ex- 
ternal field h is increased. The blue and red dots indicate the 
quantum states represented by the peak signals (cf. Fig. [T]). 



state from the rest. Following, measurement of any ob- 
servable, and even quantum state tomography, can be 
performed for the resulting ground state. With the com- 
plete knowledge of the eigenstates and the eigenvalues, we 
can in principle obtain all accessible information about 
the ground state properties; therefore, this procedure 
solves the ground-state problem when trial wave func- 
tions are available. 



V. EXPERIMENTAL PROCEDURE 

The experiments were carried out at room tempera- 
ture on a Bruker AV-400 spectrometer. The sample we 
used is the -^^C-labeled Diethyl- fluoromalonate dissolved 
in ^H-labeled chloroform. This system is a three-qubit 
quantum simulator using the nuclear spins of ^^C and 
"'^H as the system qubits to simulate the Heisenberg spins, 
and the ^^F as the probe qubit in the phase estimation al- 
gorithm (see Fig. [l]3). The internal Hamiltonian Hnmr 
of this system can be described by the following: 



H 



NMR 



= E 

jG{a,6,c} 



(6) 



where Vj is the resonance frequency of the jth spin and 
Jjk is the scalar coupling strength between spins j and 
k, with Jah = 160.7 Hz, = -194.4 Hz, and Jac = 
47.6 Hz. The relaxation time Ti and dephasing time 
T2 for each of the three nuclear spins are tabulated in 

The experimental procedure consists of three main 
parts: 1. State initialization (preparing the system qubits 
as IV^*), probe qubit as |0)), 11. Eigenvalue measurement 
by iterative phase estimation, and HI. Quantum state 
tomography. The state initialization part is rather stan- 
dard and we leave the details of it to the Appendix [31] . 
Part II is implemented with a quantum circuit as de- 
picted in Fig. [2] (see the Appendix [31] for the detailed 
circuit construction). The probe qubit is measured at the 
end of the circuit (see also Eq. ([5|). 

The resulting Fourier spectra for various cases are 
shown in Fig. [3] The positions of the peaks indicate the 
eigenvalue of the Hamiltonian H. Although the peaks 
look sharp, the errors are in fact about 22%. However, 
we are able to reduce the errors to less than 0.003% (see 
Fig. [4| by five steps of the iterative phase estimation al- 
gorithm which is described below. 



VI. ITERATIVE PHASE ESTIMATION 
ALGORITHM (IPEA) 

To improve the resolution of the energy eigenvalues, 
the information stemming from long time evolution of the 
simulated state is needed [20 . Fortunately, the required 
resources can be significantly reduced by the IPEA ap- 
proach. This is due to the symmetry of the Hamiltonian: 
since all the terms in the Hamiltonian (Eq. ([3|) commute 
with each other, they can be simulated individually, i.e.. 



iHt 



y e 



iJI^I 



(7) 



for all times t. The last term e~*^(^^^^^)^ corresponds 
to two separate local rotations, whose implementation 
is straight-forward (see the Appendix [31 ). The other 
terms e~*^^«^«^ are equivalent up to some local unitary 
rotations, and their eigenvalue spectra of which 
are 4 and —4, are the same; the eigenvalues are sym- 
metrical about zero. This means that, in order to simu- 
late each term for a time interval t, we can always find 



a shorter time r such that e" 



where 



j<feG{a,6,c} 



t = Snir/ J + r for some non- negative integer n which is 
determined by the condition: < Jr < Stt. 

Now, denote the eigenvalue, ujk = 27rJx O.X1X2X3..., by 
a string of decimal digits {xi, X2, X3...}. The first digit xi 
can be determined by a short time evolution by a probe 
qubit described in Eq. ([5|. Once xi is known, the second 
digit X2 can be iteratively determined by simulating the 
evolution for ten times longer than the previous ones: 

10 X ujkt = 2iTJt X xi.O + 2iTJt X O.X2X3... (8) 
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FIG. 4: (Color online) Experimental results of the iterative phase estimation algorithm (IPEA) for improving the accuracy of 
the measured ground- state energy, (a) and (b) There are five iterations performed; each of them improves one digit of accuracy 
in the eigenvalues. For example, consider the ground state, after the first iteration (red curve), the peak lies between —0.1 and 
—0.2; this means that the first digit of eigenvalue should be —0.1. After five iterations, the value of ground-state energy is 
determined to be —0.11936(3), with a precision of 10~^ in units of 27rJ. (c) A table listing the improvement of the numerical 
values (digits in red represent uncertainty), (d) Graphical visualization of the results in (c). The theoretical curve results from 
the improvement of the precision by a factor of 1/10 for each iteration. 



Note that the first term on the right hand side is known. 
The second term is now amplified, and can be resolved by 
the probe qubit. This means that the eigenvalue Uk can 
then be determined to two digits of precision. By repeat- 
ing this scheme iteratively for xs and so on, the eigen- 
value ujk can be determined subsequently for one digit 
after the other (cf. Fig. [4|. The accuracy of the eigen- 
values is improved from about 22% to about 0.003%. We 
note that in the IPEA performed in Refs. [23, 32 , the 
final unitary matrices are decomposed directly for each 
value of t. Therefore, one can in principle determine the 
eigenvalues to any arbitrary accuracy. However, the re- 
sources required for decomposing the unitary matrices 
grow exponentially with the system size; the methods 
implemented there are certainly unrealistic for larger sys- 
tems. Here, we exploited the symmetry of the Hamilto- 
nian, and simulate the time evolution without performing 
the decomposition of the unitary matrices. The accuracy 
of the IPEA is limited by some natural constraints. The 
details about the limitation of this method are discussed 
in the Appendix [3T| . 



VII. RESULTS AND DISCUSSION 

Once the two eigenvalues {Eq and Ei) are accurately 
determined by the IPEA, we can identify the eigenvectors 
(ground state |eo) and excited state |ei)) by the same 
quantum circuit as shown in Fig. [2^. The difference is 
that, the time r, in the controlled rotation U (r) = e"*^'^ 
is chosen to be r = tt/ (^i — Eq). This allows us to 



obtain the following state, 

-^(|eo)|0)-|ei)|l)) . (9) 

This state is very similar to the one discussed in Eq. ([2|. 
The important point is that, now each eigenstate is 
tagged by the two orthogonal states of the ancilla qubit, 
and can be determined separately, e.g. through quantum 
state tomography. 

To obtain the state in Eq (|9|, starting from the prod- 
uct state 1?/^*) |0), we first prepared the probe state as 
a superposition state with a phase e'^^^^ "preloaded" in 
it, i.e., (|0) +e^^o^ |1)) /V^. Next, after applying the 
controlled- /7(r) to the trial state = ao |eo) + ai |ei), 
we have, 

-^ao |eo) (|0) + |1)) + -^ai \e,) (|0) + e'^ |1)) . (10) 

Subsequently, we apply a single-qubit rotation gate 
i?^(-7r/2), which maps (|0) + |1))/V^ |0) and 

(|0) — |iy) — we then obtain the final state 

in Eq. 

Finally, the standard procedure of quantum state to- 
mography [34] was performed on the final states (Eq. ([9|) 
for the cases h = 0, h = 0.75/ic, and h = 1.25/ic, shown 
respectively in Fig. [s] (b)-(d). The corresponding results 
of the ground state (i.e. the |eo) part in Eq. ([9|) are 
shown in Fig. [s] (e)-(g). These density matrices allow us 
to obtain all information about the experimentally de- 
termined ground states. Fig. [6^ shows the improvement 
of the magnetization M of the final states, as compared 
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FIG. 5: (Color online) Experimental results from the quan- 
tum process tomography procedure (real parts are shown, 
imaginary parts shown in the Appendix [Slj)- (a) the ini- 
tial state IV^*). (b),(c), and (d) Three final states (Eq. ^) 
for the cases, respectively, h — h — 0.75hc, and h = 1.25/ic- 
(e),(f), and (g) The first 4x4 section of each density ma- 
trix above (after re- normalization), in the subspace where the 
probe qubit is projected to the |0) state. 



with the initial state. The inset figure shows that the 
magnitude of the deviations (blue bars) from the theo- 
retical values are always smaller then that (red bars) of 
the trial state. 

The quality of the final state pexp in the experiment is 
quantified by the fidelity F = (eo| pexp l^o 

) (cf. Eq. ffl), 

and the projection [35] P = Fj^/Q^ where Q = Tr(p^) 
is the purity of pexp- The results are shown in Fig. l6p. 
Note that the reduced density matrices (e),(f),(g) have 
better fidelities than that of the original density matri- 
ces (b),(c),(d). In Fig. [7| the weights (probabilities) of 
the eigenstates of R in the final states are shown. Note 
that, as mentioned above, the trial state captures only 
two eigenstates. Due to experimental errors, other eigen- 
states also showed up in the spectral decomposition. This 
contributes to the deviation of the magnetization (M = 
for the singlet state) as well. Note that the singlet-triplet 
switching (cf. Fig. [T]), i.e., from Fig. [t]^ to[7]i, is reliably 
captured. 

In this experiment, we are able to determine the eigen- 
values to a very high accuracy, using the iterative phase 
estimation algorithm (IPEA). The major source of errors 
(about 10% of the fidelity) of the experiment comes from 
the second step of the procedure where the overall pulse 
sequence to construct the final state Eq. ([9| is lengthy, 
and therefore is dominant ly a error; the time spent for 
this operation is about 1/10 of T2 (see the Appendix [31 ). 
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FIG. 6: (Color online) Results from the quantum state tomog- 
raphy, (a) Magnetization + (J^) for the initial states (red 
dotted line) and the final states (crossed circles) for /i = 0, 
h — 0.75hc, and h = 1.25/ic- The inset shows the absolute er- 
rors of the initial state (red bars) , and the three experimental 
values (blue bars), (b) The ground state fidelity (green) and 
projection (yellow) for the experimentally determined states 
(a)-(g) in Fig. [s] The fidelity of the initial state (blue) is 
included for comparison. 



Additionally, other errors come from the measurement 
(tomography) errors, and the inhomogeneity in the RF 
pulses and the external magnetic field. If these factors 
can be overcome, a further increase of fidelity is possible 
by using the final state of this experiment as the input 
state for another iteration of the similar distillation pro- 
cedure (see the Appendix [31 for details). 



VIII. CONCLUSION 

We have experimentally demonstrated a method to 
solve the quantum ground-state problem using an NMR 
setup. This is achieved by distilling the exact ground 
state from an input state, which has 50% overlap with 
the ground state. The eigenvalues were determined to a 
precision of the 10~^ decimal digit, after five iterations 
of the phase estimation procedure. Then, the final states 
are distilled to high values of fidelity. The method we 
developed in this experiment is scalable to more general 
Hamiltonians, and not limited to NMR systems. This 
result confirms that variational methods developed for 
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classical computing could be a good starting point for 
quantum computers, opening more possibilities for the 
purposes of quantum computation and simulation. 
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Appendix A: State initialization 

In this experiment, we used a sample of the ^^C-labeled 
Diethyl- fluoromalonate dissolved in the ^H-labeled chlo- 
roform as a three-qubit computer, where the nuclear 
spins of the ^^C and the ^H were used as the system 
qubits, and that of the ^^F was used as the probe qubit. 



U(0.195n) 



Qubit b 



n -Ti/2 



Qubit a 



■ If' 



n/2 0.195TT 

I 



Y- pulse X-pulse U(e)=exp(-ie IJJ 

I I i 

FIG. 8: (Color online) Pulse sequence for generating the input 
state IV^*). 



The structure of the molecule is shown in Fig. [T^ of the 
main text, and the physical properties are listed in the 
table of Fig. [1]d. 

Starting from the thermal equilibrium state, we first 
created the pseudo-pure state (PPS) 



Pooo = (l-e)V8 + e|000) (000| 



(Al) 



using the standard spatial average technique. Here, 
e~10~^ quantifies the strength of the polarization of the 
system, and I is the 8x8 identity matrix. Next, we 
prepared the probe qubit to the state ^(|0) + by a 
pseudo-Hadamard gate i?^(7r/2), where. 



Ri {0) = e-'^''- 



(A2) 



Here, a = x^y^ is a rotation operation applied to the 
qubit j. 

Finally, the system qubits are prepared to the initial 

state, 

1^.) = ^(101) -110)) + -^ 111) , (A3) 

by applying two single-qubit rotations and one 
controlled-rotation. The pulse sequence employed fol- 
lows: 

^ ^ i?°(0.1957r), (A4) 

where the unitary evolution, 

W''{0) = e-'^'i'" (A5) 
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is generated from the natural evolution between qubit j 
and k. 



Appendix B: Quantum circuit diagram for 
simulating the controlled-[/(t) 



Qubit c 
Qubit b 
Qubit a 



W 



rsi 

■M 
> 



rsi 

■M 
> 



FIG. 9: (Color online) Quantum circuit diagram for simulat- 
ing controlled- and controlled- 1/^^. 



(6 
c 



^^^^^^^^^^^^^^^^^^ 


t=0 

t— 1 . D/ J 




r 

t=6.4/J 







\ I I I r 

8800 8600 8400 8200 Hz 

NMR Frequency 

FIG. 11: The signals of the experimental spectra for the case 
/i = 0, where t = 0, 1.6/ J and 6.4/J. 



(a) U(ti) 
Qubit c 



-Tl/2 



Qubit b - 
(b) 

Qubit b 
Qubit a 



Tl/2 Tl/2 



-Tl/2 n/2 
U(t/2) ■ 



-Tl/2 



n/2 



Y- pulse X-pulse U(0)=exp(-ie IJJ 

. I i 

FIG. 10: (Color online) Pulse sequences for simulating 
controlled-e"^^^' and Kc(t/2). 

The controlled-^/ {t) in the phase estimation algorithm 
(see Fig. [2|l) is implemented in the following way: since 
all the terms in the the Heisenberg Hamiltonian, 



H = j [i^jl + ipl + itil) + h (/» + 4^ 



(Bl) 



commute with each other, we decompose the time evolu- 



tion operator T{t) 



-iHt • 



into three parts: 



T{t)=VAt)Vy, {t)L, (t) 



where 



L. it) 



(B2) 

(B3) 
(B4) 
(B5) 



The quantum circuit diagram for simulating the opera- 
tions controlled- and controlled- is shown in Fig.|9] 



To simulate controlled-ya;(t), we set, 

V{t/2) = 14(^2) and Wy = e"^^^^ . (B6) 
(alternatively, i^); to simulate controlled-l^2(^)7 we set 

V{t/2) = Vy,{t/2) and W^, = e''^'^ . (B7) 

Note that the control is "on" when the probe qubit is 
in the |0) state. In this case, the first three quantum 
gates cancel the last gate V(t/2), making it effectively 
an identity gate. When the controlling qubit is in the 
"off' state, this circuit executes two V{t/2) gates. 

The pulse sequences for generating the controlled- 
e-*^^x gate are, 

and the pulse sequences of the Vx{t/2) gate is: 

The corresponding diagrams of the pulse sequence are 
shown in Fig. [lOj 

Appendix C: Measurement of the probe qubit 

Here we explain the measurement method of the NMR 
signal of the probe qubit (see Eq. ([5|). Denote the off- 
diagonal elements of Pprobe{t) as. 



(CI) 



The phase shift (j)t can be obtained by using the method 
of quadrature detection which serves as a phase detector. 
By measuring the integrate value of the peak in NMR 
spectrum, we can obtain the value of \Mt\. 
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To calibrate the system, we adjust the phase of the 
NMR spectrum such that becomes the reference 
phase, and normahze its peak intensity as 1. Some of the 



experimental data of the spectra are shown in Fig. 11 for 
the case of = 0, at t = 0, 1.6/ J and 6. 4/ J. 

By simulating the Hamiltonian evolution for different 
times, a range of frequency spectrum of \Mt \ e*"^* can be 
obtained by the method of discrete Fourier transforma- 
tion (DFT). The Fourier-transformed spectra are shown 
in Fig. |3] for the cases oi h = 0,0.75/ic7 and 1.25hc^ re- 
spectively. For each spectrum, totally 128 data points 
were collected. 



Appendix D: The precision limit of the iterative 
phase estimation algorithm 

In principle, it is possible to simulate the time evo- 
lution for an arbitrarily long time by mapping it back 
to a corresponding short time evolution. In practice, 
this method is limited by the precision of J, which is 
determined independently in the experiment. Here we 
show that when J is changed by a small amount, i.e., 
J ^ J -\- 5J ^ then the error for determining the phase 
angle for the mapping goes as 8n7r x (SJ/J). In this ex- 
periment, we are able to determine the eigenvalues to the 
fifth digit of accuracy (see Fig. |4|. 

To elaborate more, let us consider the iterative phase 
measurement. For the moment, let use consider one of 
the terms in the Heisenberg Hamiltonian H de- 

fined in Eq. ([3|. We want to find the value of a such 
that. 



(Dl) 



where t = Snyr/J + r, and n is determined by the condi- 
tion that < Jr < Stt. Ideally, we have, 

a = Jt = J X {Snn/J) ^ Jr = Jr (ideal). (D2) 

If there is a fluctuation of J ^ J + (5 J, then the wrong 
a , call it as , is: 



as = {J ^ SJ) t = Jr^SJ X {Snn/J + r) 



(D3) 



The change of the phase angle, Aa = a — as = SJ x t/is 
therefore equal to 



Aa = SJ X (8n7r/J + r) 



(D4) 



which becomes Aa ^ Snir x (SJ/J) for large n. 

In the phase estimation algorithm for determining the 
eigenvalue if we set Et ^ a (up to some constant), 
then SE ^ Aa/t ~ SJ. In this experiment, SJ/J ~ 
0.01%, which makes SE/E ^ 0.01%. This is in agree- 
ment with the data of the ground-state and excited-state 
energies m Fig. m the accuracy is about 10 ^ x 27r J. 

On the other hand, we comment one point which may 
need attention in the implementation of the iterative 



phase estimation algorithm described in this work. In 
our method, although we have an accuracy about 0.04 (in 
units of 27rJ) in reading the digit for every iteration, it 
is not guaranteed that the digit determined is correct for 
all iterations; some error-correction procedure is needed. 
This is because in some exceptional cases, for example, in 
the second iteration of the excited state, the experiment 
result is 0.038916 and the theoretic value is 0.039788. 
If, unfortunately, we obtained the experiment result as 
0.041788 instead, which has a difference of about 0.002 
from theoretic value, in our procedure, we would con- 
clude that the second digit of the energy is 4, but the 
right answer is 3. We can only solve this problem in the 
following iteration; in the next iteration, even if we used 
the wrong value of the second digit, we will obtain a peak 
not lying between and 1. So we can determine that the 
second digit should be 3 instead. 



Appendix E: Generalization to the cases of multiple 
eigenvalues 

In this experiment, we have chosen the case of the trial 
state lip^) that captures two out of four eigenstates of 
the two-spin Hamiltonian. Therefore, we can use a single 
qubit (two states) to resolve the two distinct eigenvalues, 
and map the final state into the form defined in Eq. ([9|, 
which is then analyzed by a quantum state tomography 
to extract the information about the ground state |eo). 

In general, a trial state may capture more than two 
eigenvalues. In this case, our procedure needs to be gen- 
eralized. However, there is nothing fundamentally new, 
except for a more laborious repetition of the same pro- 
cedures. This is the reason we decided to work on the 
specific case of the trial state being the linear combina- 
tion of two eigenstates only. 

To explain the details of how it works, we assume the 
ground-state energy of H is unique. Define the first ex- 
cited state as |ei). Then, any trial state can be decom- 
posed into the following form: 



IV^*) = tto |eo) + cii |ei) + a2 \e2) 



(El) 



where |ao| + + |a2| =1, and |e2) represents the 
linear combination of all higher energy states captured by 
IV^*). Then, we perform the phase estimation algorithm, 
using a single probe qubit (cf. Eq. ([5|), and obtain all 
of the eigenvalues. Performing the same procedure for 
getting Eq. ([9|, we can obtain the following state: 

{bo \eo) + b20 |e2» |0) + (61 |ei) + 621 |e2)) |1> , (E2) 

where 1 60 1 ^ + 1 61 1 ^ + 1 620 1 ^ + 1 ^21 1 ^ = 1 • Now, if we perform 
a state tomography, and extract the first part of the state, 
we obtain a new state 



^0 |eo) + ^20 |e2) 



(E3) 



which contains no eigenstate |ei). If we use this new 
state as the new trial state for another cycle, we get one 
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(a) Input state (b) h = (e) h = 




less eigen-energy to worry about. Therefore, we can in 
principle eliminate the higher eigenstates one after one, 
and obtain the ground state in the end, using a single 
probe qubit. 



FIG. 12: (Color online) Imaginary parts of the tomography 
results, (a) the initial state IV^*). (b),(c), and (d) Three final 
states (Eq. Q) for the cases, respectively, h — 0, h = 0.75/ic, 
and h = 1.25/ic. (e), (f), and (g) The first 4x4 part of each 
density matrix above. 



Experimental procedures Running time 

State initialization ~ 10 ms 

Fourier spectrum for the case h=0 ~ 120 ms 

Fourier spectra for cases h?iO ~ 150 ms 

Distillation + State tomography ~ 100 ms 



Appendix F: Supplementary data 



1. In Fig. [T2j the imaginary parts of the results from 
the quantum state tomography are shown. 



2. In Fig.[T3j the running times of various experimen- 
tal procedures are shown. 



FIG. 13: (Color online) Running times for various experimen- 
tal procedures. The time scale is taken from the longest ones. 
Distillation refers to the procedure to obtain Eq. ([9|. 



